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ABSTRACT 

High cycle fatigue of metals typically occurs through long term exposure to time varying 
loads which, although modest in amplitude, give rise to microscopic cracks that can 
ultimately propagate to failure. The fatigue life of a component is primarily dependent on the 
stress amplitude response at critical failure locations. For most vibration tests, it is common 
to assume a Gaussian distribution of both the input acceleration and stress response. In real 
life, however, it is common to experience non-Gaussian acceleration input, and this can cause 
the response to be non-Gaussian. Examples of non-Gaussian loads include road irregularities 
such as potholes in the automotive world or turbulent boundary layer pressure fluctuations for 
the aerospace sector or more generally wind, wave or high amplitude acoustic loads. 

The paper first reviews some of the methods used to generate non-Gaussian excitation 
signals with a given power spectral density and kurtosis. The kurtosis of the response is 
examined once the signal is passed through a linear time invariant system. Finally an 
algorithm is presented that determines the output kurtosis based upon the input kurtosis, the 
input power spectral density and the frequency response function of the system. The 
algorithm is validated using numerical simulations. Direct applications of these results 
include improved fatigue life estimations and a method to accelerate shaker tests by 
generating high kurtosis, non-Gaussian drive signals. 
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1. INTRODUCTION 


Any operational loads monitoring exercise will show that loads differ from one application to 
another. Some loadings, like shocks, are of a transient nature, while others, like turbulent 
boundary layer, are of a random nature. A random process can be stationary or non- 
stationary, Gaussian or non-Gaussian. A process is referred to as ‘stationary’ if its probability 
distribution does not change with respect to time. A random variable x is said to be Gaussian 
if its probability density function (PDF) is given by 
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where ju is the mean and a is the standard deviation. 

A non-stationary process is often represented as a signal with time-varying variance. Gust 
loading is an important part of the fatigue induced on an aircraft structure. In this case, gusts 
with various amplitudes and directions must be considered. Similarly, road irregularities and 
vehicle speed variations are the primary cause of failure in automotive structures and 
components. Both gust loading and road irregularities can be represented by signals with 
time -varying variance. 

It is desirable for a design engineer or a test technician to have the means to compute the 
kurtosis of the stress responses at critical locations of a component or structure. This would 
allow the design engineer to perform more accurate fatigue life estimations using a spectral 
fatigue approach and the test technician to know the acceleration factor that can be achieved 
by using the existing kurtosis control capabilities of a shaker. 

Kurtosis is one of the principal metrics describing non-Gaussian features of a unimodal 
probability distribution. It is often described as a measure of the "peakedness" of a process . 
The kurtosis is calculated from the fourth central moment (m 4 ) of a random variable x, 
normalized by the square of the variance (w 2 ), i.e., 
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Kurtosis of a Gaussian distribution is 3.0. A process is said to be leptokurtic if its kurtosis is 
higher than 3.0, and platykurtic if smaller than 3.0. 

In recent work by Rizzi et al. [1], it was found that a leptokurtic loading having a steady 
rate of instantaneous, high-excursion peaks produced essentially the same response as if the 
load was Gaussian. In contrast, the response to a non-Gaussian loading having the same 
kurtosis, but with bursts of high values was found to be non-Gaussian and also leptokurtic. 
Over a practical range of damping, it was found that the linear response to a non-Gaussian 
loading was Gaussian when the period of the system impulse response is much greater than 
the rate of the peaks in the loading. A lower damping would typically reduce the excess 
kurtosis in the response. 

In this paper the previous study is further extended by concentrating more on stresses and 
fatigue damage. Also, the processes for simulating a non-Gaussian signal are more fully 
discussed. A formulation is proposed that gives the statistics of the response due to a steady 
non-Gaussian loading. Finally, the case of non-Gaussian loading with bursts of high 
amplitude excursions is further explored and the relationship between the system impulse 
response and the frequency at which the bursts occur in the loading is illustrated. 


2. SIGNAL GENERATION 

In the following sections, three forms of loading signal generation are undertaken: Gaussian, a 
non-Gaussian form using a PDF transform approach, and a non-Gaussian form using 
amplitude modulation. 


2.1 Gaussian stationary signal generation 

Time histories corresponding to particular power spectral densities (PSD) may be generated in 
various manners. One method is based on time domain filtering of white noise. In this case, 
the coefficients of the finite impulse response filter are calculated from an amplitude specified 
by the square root of the PSD, and a phase specified by a uniform random distribution 
between [-7r, 7r]. 

Another method of signal generation is based on the inverse Fourier transform (IFT) of a 
frequency domain signal constructed with an amplitude specified by the PSD combined with 
a random phase. Newland [2] expresses the time signal reconstruction of a stationary, zero- 
mean process as: 

x{t) = Re(^2AfG xx (kAf)e 2MAf e ^ ) (3) 

k= 0 

with G xx the PSD to be inverted, (p the randomized phase between [- 71 , n], A f the frequency 
resolution and N the number of frequency lines. 

The length and sample rate of the signal generated depend on the frequency resolution and 
the maximum frequency of the PSD respectively. Note that the signal obtained with the IFT 
method is generally considered as having a lower level of randomness than filtered white 
noise [3]. For instance a filtered white noise is expected to be unbounded, whereas the 
maximum of a signal obtained with the IFT method is limited as shown by: 

max(x(t))<|^/2AfG xx (kAf) (4) 

k=0 

An extension of the previous method for generating a colored Gaussian process is the so- 
called ‘inverse Welch’ method, as it was modeled after the Welch method for estimating 
spectral density [3]. Frames of random data obtained from the IFT method are multiplied by a 
window, overlapped, and added to form a nearly stationary output. The weighting windows 
allow continuous frequency content and help limit leakage. A typical window is the Hanning 
window. The level of overlap dictates whether the resulting signal will be stationary or not. 
For instance, Hanning windows, which are essentially cosine functions, when overlapped by 
67% will provide a uniform statistical weighting which in turn will generate a stationary 
signal. Figure 1 shows a target PSD overlaid with the PSD of the signals obtained after 
filtering white noise and using the inverse Welch method. 



Figure 1 : Comparison of prescribed PSD (red) with PSD of generated signals using the white 
noise filtering (blue) and inverse Welch (green) approaches. 

Both the filtered white noise and the inverse Welch approaches give very similar and 
satisfactory results. The inverse Welch method is adopted here because it offers more 
flexibility, particularly when non-stationary signal generation is required. 


2.2 Steady non-Gaussian signal generation - PDF transform technique 


The technique introduced in this section can be used to obtain a non-Gaussian stationary 
signal, showing a steady rate of instantaneous peaks. 

The PDF of a random variable x can be transformed to the PDF of a different random 
variable z using the PDF transform [4, 5]: 
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where there is a one-to-one monotonic relationship between x and z, that is, z = g(x). Several 
formulae for the transformation g(.) have been proposed in the literature [ 6 , 7]. The idea is 
based on the fact that squaring part of a Gaussian signal will have the effect of skewing it, 
whereas cubing a Gaussian signal will make it leptokurtic. Several authors [ 6 , 7] thought of 
using a polynomial transform of the general form: 

Z = c + aX + /3X 2 + yX i ( 6 ) 


with c, a, ft and y coefficients driving the type and amount of non-Gaussianity that is required. 
Typically, the constant c in equation ( 6 ) specifies the mean value, while a drives the standard 
deviation of Z. The ft and y coefficients influence the skewness and kurtosis respectively. 
The goal is therefore to find a relationship between the coefficients in equation ( 6 ) and the 
required statistics for Z. 

Winterstein [ 6 ] considered using Hermite polynomials to expand a Gaussian PDF to match 
the first statistical moments required of the transformed signal z. For a leptokurtic process, 
Winterstein’ s model [ 6 , 8 ] defines the direct transformation g that allows one to obtain z(t) 
from x(t) as: 
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where y 3 and 74 are the expected skewness and expected kurtosis respectively. The 
Winterstein model was adopted here for the generation of steady non-Gaussian signals 
because it appears to be widely accepted and it is known to provide an accurate representation 
of a wide range of nonlinear behaviours [ 8 ]. Note that this method has one major drawback: 
it introduces harmonic distortions which reduce the dynamic range in the spectral density of 
the generated data [9]. 


2.3 Burst non-Gaussian signal generation - amplitude modulation technique 

Another class of high kurtosis random signals is made of non-stationary processes. Although 
vibration is generally of random form, it is often non-stationary. In the automotive world, this 
arises from changes in road surface quality and vehicle speed. A non-stationary process can 
be represented as a signal with time-varying variance, typically obtained through a 
modulation process. 

The process for generating a non-Gaussian signal with bursts starts with the generation of a 
stationary Gaussian signal. This signal is then amplitude modulated by a low frequency 
carrier wave, independent of the Gaussian signal. In order to create the large excursions, the 
amplitude of each cycle of the carrier wave is a random variable. As described by Smallwood 
in [3], aft distribution is chosen for the realizations of the random variable. The ft distribution 
is chosen because it generates only positive values and is extremely flexible. The mean value 
of the distribution is fixed to 1.0, while its variance is used to control the kurtosis. Note that 
the skewness obtained using this method is zero. 


The non-Gaussian signal z(t) can therefore be written z(t) = a(t).x(t), where a(t) is a 
amplitude modulation coefficient that changes over time and x(t) are realizations of a 
stationary Gaussian process. The bursts that are created in the resulting signal occur at a 
given rate, which derives from the period of the modulating signal a(t). Also, the statistics of 
z(t) can be derived from the statistics of the independent signals a(t) and x(t). 

It is therefore possible to construct a burst non-Gaussian loading with a prescribed PSD 
and kurtosis value by generating two signals that are multiplied. First, a Gaussian stationary 
signal is generated from the given PSD. Secondly, a low frequency varying amplitude carrier 
wave is constructed. The amplitudes of the waves generated in this study are realizations of a 
P distribution with parameters carefully chosen such that the kurtosis of the wave is one third 
of the required kurtosis. Finally, the resulting signal is obtained by multiplying both signals 
and then scaling to the original RMS value. Such a burst non-Gaussian signal is illustrated in 
Figure 2. 



(1) Gaussian 


(2) Carrier wave 


(l)x(2) 

Burst non-Gaussian 


Figure 2: Construction of a burst non-Gaussian signal 


3. NUMERICAL SIMULATIONS: LINEAR TIME-INVARIANT SYSTEM 

The response kurtosis, maximum stress and fatigue damage for a simple, but realistic, 
engineering structure is next considered to demonstrate how these quantities change for non- 
Gaussian loadings with the same kurtosis but with differing temporal character. 


3.1 The Finite Element model, boundary conditions and material properties 

The engineering structure considered is an off-the-shelf automotive horn bracket. The 
finite element model of the bracket is shown in Figure 3. The bracket was modelled with the 
finite element analysis program ABAQUS [10]. The bracket is constructed with 168 S3R and 
5382 S4R shell elements. Two lumped masses representing the horns were constructed with 
8074 C3D4 constant stress tetrahedron solid elements. Different mass values can simulate 
asymmetry of the weight distribution. The bracket was constructed from carbon steel, with 
the properties indicated in Table 1. 

The fatigue S-N curve, when plotted on a log-log scale, is specified by a line with constant 
slope, and no fatigue limit. Under simulated test, the specimen is clamped to the test rig at the 
left end and excited in the transverse direction, as shown in Figure 3. The following section 
describes the applied loading. 


3.2 Input signal generation 

The prescribed acceleration PSD, shown in Figure 4, originates from MIL-STD-810F [11]. 
The loading has a frequency range of 10-2000 Hz and an amplitude of 3.7 g RMS. 





Figure 3: The finite element model of the bracket with boundary conditions. 



Property 

Value 

Mechanical properties 

Young Modulus 

210 GPa 

Density 

7890 kg/m J 

Modal damping 

0.025 

Yield Stress 

400 MPa 

Ultimate tensile strength 

700 MPa 

Fatigue properties 

Basquin slope 

14 

Stress Range Intercept 

1200 MPa 


Table 1: Material properties of hom bracket. 



Frequency (Hz) 

Figure 4: The excitation PSD used for Gaussian and non-Gaussian loadings. 

Three different excitation signals were produced from this PSD. The first one is a 
stationary Gaussian signal, whereas the other two are leptokurtic with a high kurtosis value of 
roughly 7. All signals share the same PSD and therefore the same RMS value. The first 
leptokurtic signal was obtained by the PDF transform as described in section 2.2, while the 
latter results from the amplitude modulation described in section 2.3. For the latter, the 
modulator is made of a series of 50% overlapped Hanning windows, each comprised of 1024 
points. All signals are 50 seconds in duration and are sampled at 4096 samples per second. 
The signals and their PDFs are illustrated in the Figure 5 and Figure 6. In Figure 5, the steady 
rate of high amplitude peaks is shown in the non-Gaussian signal generated using the PDF 
transformed technique, whereas the signal generated with the amplitude modulated approach 
shows bursts. The same kurtosis for both non-Gaussian signals is clearly seen in the PDFs in 
Figure 6. 
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Figure 5: The excitation signals generated: Gaussian (red), steady non-Gaussian using PDF 
transform method (blue) and burst non-Gaussian using amplitude modulation method (green). 



Figure 6: The PDFs of the excitation signals: Gaussian (red), steady non-Gaussian (blue) and 

burst non-Gaussian (green). 


3.3 Dynamic finite element analyses 

Two linear dynamic response analyses were performed via the finite element method: a forced 
frequency response analysis and a transient analysis. Results from the forced frequency 
response analysis serve as input to a spectral fatigue calculation, while results from the 
transient analysis serve as input to a fatigue analysis performed using time domain data and a 
rainflow cycle [12] counting scheme. 

The forced frequency response analysis was performed by exciting the structure in physical 
degrees-of-freedom with a 1 g vertical sinusoidal base acceleration input over the frequency 
range of interest. The resulting frequency response function (FRF) between the loading and 
the response is obtained in this manner. An example FRF for the absolute maximum principal 
stress (Pa/g) is shown in Figure 7 for the critical node 3100 (see Figure 10). The absolute 
maximum principal stress is defined as the principal stress with the largest magnitude. Once 
the FRF is known, the response PSD Gyy can be calculated for any loading PSD Gxx as: 

G yy (m) = |/7(m)| 2 -G zy (m) (8) 

where H(co) represents the FRF. For the fatigue analysis, stress cycles are derived statistically 
from the stress response PSD, as described in the next section. 




Frequency (Hz) 

Figure 7: Stress FRF obtained at a critical node (3 100). 

For the time domain approach, a modal transient analysis is performed in ABAQUS to 
determine the modal response time histories for each loading condition. Stress time histories 
were obtained by modal superposition of the modal stresses with the modal coordinates using 
the commercial software nCode DesignLife™. This approach allows the size of the finite 
element result file to be kept small, since stress time histories at all nodes are not required, 
and those of interest can be generated quickly. The accuracy of the modal approach is 
dependent on the number of modes included in the analysis. In this study, the first 8 modes 
were used, covering natural frequencies up to 300 Hz. This selection was made based on the 
FRF at the critical node, shown in Figure 7, where it is seen that the main excited modes are 
all within the 0-300 Hz range. 

3.4 Fatigue analyses 

Fatigue damage can be calculated using frequency or time domain data [13, 14]. Both 
approaches are considered here in order to highlight differences between them. 

Frequency Domain Fatigue Analysis 

A common approach for characterizing random vibrations is to compute the ensemble 
averaged PSDs. This representation is useful because it enables identification of dominant 
frequency components, as well as the average energy level within each frequency band. With 
this approach, however, the PDF of stress cycle rainflow ranges must be derived from the 
stress PSD in order to estimate the fatigue life. There are a number of algorithms for doing 
so, several of which are supported by nCode DesignLife used in this work. The particular 
method used was given by Lalanne [15] and is based on the seminal work of Rice [16]. It is 
the most often used and can handle a wide range of loading conditions [13]. Rice 
demonstrated that for a signal of arbitrary bandwidth, the PDF of the peaks could be obtained 
from the weighted sum of the Rayleigh and Gaussian distributions. Lalanne reasoned that 
over a sufficiently long period of time, the PDF of rainflow ranges would tend to the PDF of 
the peaks and thereby demonstrated that Rice’s original formula would also suffice for 
rainflow ranges. 

The spectral fatigue method is robust, efficient and offers advantages over the time domain 
approach (next described) when limited stress time history data is available. However, this 
approach is only valid when the underlying random process is ergodic, stationary and 
Gaussian [13]. Here, an ergodic process is meant to describe a random process in which the 
time average of one sequence of events is the same as the ensemble average (and vice versa). 

Time Domain Fatigue Analysis 

In the time domain, time histories of stresses are directly counted using a rainflow cycle 
counting technique [17]. The damage from each cycle is then linearly summed using the 


Palmgren-Miner [18, 19] damage accumulation law. The time domain approach is applicable 
to any type of signal, whether it be stochastic or deterministic. This approach, however, 
proves computationally intensive for stochastic loadings because long stress time histories are 
required to generate the tails of the stress range histogram in a statistically accurate manner. 
Poorly realized extremes can have a detrimental effect on the fatigue life estimate as the most 
damaging events are attributable to the high stress ranges in the tails. Convergence in damage 
estimates therefore improve with an increasing number of data points in a signal [13]. In 
many practical cases, it is simply not possible to obtain the required duration of time signal 
for accurate rainflow counting. In these cases especially, a frequency domain approach offers 
significantly more accurate fatigue life estimates under the aforementioned conditions. 

3.5 Interpretation of the fatigue results 

Contour plots of the fatigue damage obtained using the time domain method are shown in 
Figure 8 for the three loadings described in section 3.2. The critical node (3100) is the same 
in the three cases; it is in the top right radius and is represented by a small white marker on 
the plots in Figure 8. Note that a second critical point in the top left radius exists due to 
symmetry. Under the linear cumulative damage law, failure occurs when the damage equals 
unity. Thus, an accumulated damage of 10" 8 for a 50 second exposure corresponds to an 
estimated fatigue life of 5xl0 9 seconds or roughly 158 years for that loading condition. 
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Figure 8: Contour plot of fatigue damage from various types of excitations. 


Table 2 provides further evidence that the non-Gaussian load made of bursts is indeed the 
most severe. It is clear that leptokurtic input excitation signals are more damaging than 
Gaussian signals, even though all excitation signals share the same energy level. However, it 
is also clear that damage is not only related to the kurtosis of the input. Although the steady 
and bursts non-Gaussian signal have nearly the same kurtosis value of approximately 7, the 
damage resulting from their application is very different. The signal made of bursts clearly 
generates more damage in this example. This is due to the much higher maximum stress 
range. The steady non-Gaussian signal is approximately 1.5 times more damaging than the 
Gaussian stationary one, whereas the bursts non-Gaussian is over 9 times more severe in 
terms of fatigue. A plot of the stress range PDF is shown in Figure 9. Here, the higher stress 
ranges associated with the burst non-Gaussian loading are clearly evident. 


Signal 

Input 
RMS (g) 

Input 

kurtosis 

Output 
RMS (MPa) 

Output 

kurtosis 

Max Stress 
Range (MPa) 

Max Damage 

Normal 

3.04 

3 

47.1 

2.95 

334 

5.5xl0" 8 

Steady non- 
Gaussian 

3.04 

7.11 

48.2 

2.95 

336 

8.6xl0" 8 

Burst non- 
Gaussian 

3.03 

7.16 

47.2 

3.57 

405 

5.16xl0" 7 


Table 2: Comparison of stress response and accumulated damage 
for three loading conditions. 



Principal Stress (MPa) 

Figure 9: The PDFs of the stress response: Gaussian (red), steady non-Gaussian (blue) and 

burst non-Gaussian (green). 

One can compare the time domain analysis and spectral domain analysis results to 
understand the differences between the approaches. The contour plot of damage resulting 
from the spectral domain fatigue analysis for the Gaussian case is illustrated in Figure 10. As 
expected, the spectral domain analysis gives a higher damage of 8.1xl0" 8 compared with 
5.5xl0’ 8 for the rainflow counting time domain method because the tails of the former are 
statistically better defined. A plot of the rainflow counting histograms for the Gaussian case 
using the two approaches clearly indicates this point, as shown in Figure 1 1 . 



Figure 10: Contour plot of damage from an excitation defined as a PSD. 
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Figure 1 1 : Rainflow histograms from time domain (blue) and spectral domain (red) 

for a 50 second exposure duration. 


The correlation between the two approaches increases rapidly with signal duration as 
shown in Figure 12, where three longer realizations of Gaussian signals were considered: 100, 
200 and 500 seconds It is clear that the damage computed from the time domain and spectral 
approaches converge to within 5% when the signal reaches 500 seconds. However, 500 
seconds of simulated time for the time domain analysis is more computationally expensive 
than the spectral fatigue calculation. 
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Figure 12: Comparison of damage from spectral and time approaches with duration of signal. 


Finally, recall that the spectral fatigue approach is a robust, efficient and accurate way of 
estimating fatigue life when the underlying random process is Gaussian. This limitation 
stems from the assumption of a Gaussian process, which is used when the distribution of 
peaks is computed [2, 13, 15, 16]. An extension of the spectral fatigue approach for non- 
Gaussian processes has been studied and seems to give satisfactory results [8]. However, in 
order to apply it, it is necessary to evaluate the kurtosis of the stress response knowing the 
kurtosis of the load and the FRF. This is the subject of the next section. 


4. MODELLING INPUT-OUTPUT STATISTICS 

Expressions for the kurtosis of the stress response in terms of the kurtosis of the load and the 
frequency response of the system are next discussed for the cases of steady and burst non- 
Gaussian inputs. In the following, the structural response is assumed to be linear so that the 
structure can be modelled as a linear time-invariant (LTI) system. 



4.1 Input-output kurtosis relationship in the case of a steady non-Gaussian input 


For the steady non-Gaussian case, statistics of the output signal for a LTI system can be 
calculated using the cumulants [20]. Cumulants are alternatives to the statistical moments and 
have interesting properties such as additivity and homogeneity. Additivity means that the 
cumulants of the sum of n independent random variables are the sum of the cumulants of each 
random variable. Therefore, in the case of an independent and identically distributed (i.i.d.) 
process, cumulants are additive under convolution. Homogeneity means that if a random 
variable X is multiplied by a constant c, then the m th cumulant of this product will be the m th 
cumulant of X multiplied by c m . Therefore, if the convolution is interpreted as a weighted 
sum of shifted copies of the input signal, the 2 nd and 4 th cumulants of the response signal c V 2 
and c y 4 can be derived from the 2 nd and 4 th cumulants of the input signal c X 2 and c X 4 and from 
the coefficients h t of the impulse response function as: 


C y2 C x2 ^ Hi 

C ,4 =C , 4 IX 


(9) 


The relationship between cumulants and moments are given by the Leonov and Shiryayev 
formulae [20] for the 2 nd and 4 th order cumulants c y 2 and c y 4 as a function of the 2 
order moments m y 2 and m y 4 . In the case of a zero-mean signal, these are: 
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The output kurtosis is obtained by inserting equations (10) into equation (2) giving 
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The output kurtosis of a zero-mean process expressed in terms of the impulse response and 
input kurtosis is obtained by substituting equations (9) into equation (11) giving 

ix 
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Now consider how the behavior of the output kurtosis is changed as a function of the 
characteristics of the impulse response function (IRF). The IRF corresponding to systems 
with very low damping will oscillate longer than the IRF corresponding to a highly damped 
system. In relation to equation (12), an IRF with long decaying oscillations will maximize the 
squared terms in the denominator compared to the terms raised to the 4 th power in the 
numerator. This will lead to a lower output kurtosis for a lightly damped system than for a 
more heavily damped one. This behavior is consistent with the observation that the output 
kurtosis tends to Gaussian as the bandwidth of the filter, i.e, the amount of damping, 
decreases [1, 21]. 

Equation (12) is only valid for i.i.d. processes which have no temporal dependency and are 
spectrally flat. Equation (12) can be generalized to non-spectrally flat (colored) noise by 
considering the colored noise as a filtered version of spectrally flat noise, that is, 
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where c, are the values of the autocorrelation function of the input noise and f are the values 
of the convolution of the autocorrelation function with the IRF. 


4.2 Input-output kurtosis relationship in the case of burst non-Gaussian input 

Finding a closed-form solution for the burst non-Gaussian case is somewhat more 
complicated than the steady non-Gaussian case. The approach proposed here is one of 


deduction and has been developed with the aid of numerical simulations. The numerical 
simulations were conducted by generating 20 different burst non-Gaussian signals and 
passing them through 15 filters of various bandwidths to deduce a relationship between the 
‘kurtosis rate’ and the bandwidth. Here the ‘kurtosis rate’ represents the input-output kurtosis 
relationship and is defined for k x ^ 3 as: 

_ K y ~ 3 n/n 


Simulations have shown that the kurtosis rate fluctuates between 0 and 1, which means that 
the kurtosis of the output signal varies between 3.0 and the kurtosis of the input signal. 

The behavior of the ‘kurtosis rate’ can be explained in terms of the time required for the 
system IRF to decrease to some fraction of its initial amplitude, tf and the period of the bursts, 
T a . In this work, a fraction of 0.1 was used for tf If the IRF is very short compared to the 
period of bursts T a in the input signal, then each block is filtered independently and the output 
signal will just be the filtered stationary part of the input signal, modulated with the same 
periodic waveform a(t). In this case, the kurtosis of the output signal will equal that of the 
input signal. Numerical simulations have shown that considering tf smaller than 1/1 0 th of the 
period of bursts in the input signal results in an output kurtosis similar, if not equal, to the 
input kurtosis, that is, a ‘kurtosis rate’ of 1. Thus, tf < 0.1 T a is a general criterion for the 
input kurtosis to be reproduced in the output. 

If the IRF is very long compared to the period of bursts T a in the input signal, the 
convolution process will average out the high amplitude excursions. In this case, the output 
distribution tends towards a Gaussian distribution by virtue of the Central Limit Theorem [1, 
22]. Hence, the ‘kurtosis rate’ becomes 0. It is interesting to note that a burst non-Gaussian 
signal is similar to cyclostationary noise [23], which is formed by a stationary noise x(t), 
multiplied by a waveform a(t) of period T a = ///„, where f a is the cyclic frequency of oscillator 
a(t). In this case, the stationary process x(t) is supposed Gaussian. Authors from the 
electronics industry have shown that filtering a cyclostationary noise process with cyclic 
frequency f a through a LTI band-pass filter (with bandwidth B < f a /2 ) results in stationary 
output noise [24]. The present case of burst noise is a bit different than cyclostationary noise 
as the periodic wave is generally made of variable amplitudes overlapping weighting 
windows. Here, numerical simulations have shown that considering tf > 4T a is a robust 
criterion for the output signal to become stationary and Gaussian, no matter what the input 
kurtosis was. 

Finally, if the IRF length lies somewhere between these extremes, that is, if 0.1 T a < tf < 
4T a , then a linear-logarithmic interpolation is adopted. To summarize, the ‘kurtosis rate’ law 
is expressed as: 

0 t f /T a > 4 


K. 


Ky~ 3-0 

*-,- 3.0 


log(4)-log (t f /T) 
log (4) + 1 


0.l>t f /T a <4 


(15) 


1 t f /T a < 0.1 

Equation (15) shows that the kurtosis of the output signal tends to 3.0 as the bandwidth of the 
filter decreases. Equation (15) also illustrates that the full input kurtosis can be passed to the 
output if the input signal is modulated with a sufficiently low frequency carrier. 


4.3 Comparison with simulation 

Simulations have shown different damage values for the steady and burst non-Gaussian input 
loads, although both input signals had the same energy and kurtosis. The difference in 


damage is due to the response kurtosis being higher in the case of a burst non-Gaussian 
signal. In this section, the kurtosis input-output relationships found in the two previous 
sections will be used to retrieve the statistics obtained by FE time domain simulation. The FE 
simulation is described in section 3. 

Both relationships are based on the IRF of the system. The IRF relating output stress to 
input acceleration at the critical node is obtained by taking the inverse Fourier transform of 
the FRF shown in Figure 7. The obtained IRF is shown in Figure 13. It has the classical form 
of a decaying exponential and the decay rate is driven by the amount of damping. 



Figure 13: Impulse response function obtained from stress FRF at a critical node (3100). 

Let us look first at the steady non-Gaussian loading signal, which is presented in section 
3.2. Equation (13) can be used to predict the kurtosis of the stress response and this leads to a 
predicted value of 3 compared to a ‘measured’ value of 2.95. Here, the measured value is that 
computed from the simulated time domain response calculation. 

Now consider the burst loading signal, which is also presented in section 3.2. There are 
roughly 4 bursts per seconds, which means that the period of bursts T a is 0.25. The IRF at the 
critical node takes 0.65 seconds to reach l/10 th of its initial amplitude. From equation (15), 
the kurtosis rate is 0.11 and since the kurtosis of the excitation is 7.16, the predicted kurtosis 
of the stress response is 3.45. This is close to the measured value 3.57. 

The above results demonstrate that equations (13) and (15) can be used to predict the 
kurtosis of the stress output for a LTI system exposed to steady and burst non-Gaussian 
loadings, respectively. Together with the work of Benasciutti and Tovo [8], spectral fatigue 
life estimates based on the kurtosis of the excitation and characteristics of the LTI system are 
possible. A demonstration of this is left as an area of future work. 

5. CONCLUDING REMARKS 

The fatigue life of an automotive horn bracket was investigated to determine how changes in 
the loading kurtosis influence the results. Three base acceleration loadings having the same 
PSD and RMS level, but different kurtosis were considered: a Gaussian loading and two 
leptokurtic non-Gaussian loadings of a steady and burst nature. Using results from a 
computationally intensive time domain analysis, it was found that burst non-Gaussian loading 
produced the most damage, followed by the steady non-Gaussian loading and lastly the 
Gaussian loading. 

In the case of non-Gaussian loads, the stress response may tend to Gaussian depending on 
type of non-Gaussian input signal and on the bandwidth of the system’s frequency response. 
Non-stationary processes with bursts of energy spaced in time easily propagate their 
leptokurtic nature to the response. The more bursts are separated in time, the more input 
kurtosis will be transferred to the output response. On the contrary, steady non-Gaussian 



processes essentially produce responses that tend to Gaussian, especially for narrow band 
systems. 

An alternative spectral fatigue analysis was undertaken for the Gaussian loading condition. 
The spectral fatigue analysis is considered more attractive, when applicable, than the time 
domain approach because long time histories are needed to generate statistically accurate 
stress distributions in the latter. The two approaches were shown to converge when the time 
domain simulation was carried out over an extended period. However, the standard spectral 
fatigue analysis is considered valid only for Gaussian processes. In the case of leptokurtic 
processes, it is anticipated that the standard spectral fatigue analysis under-estimates damage. 

An extension of the spectral fatigue analysis for non-Gaussian processes has been studied 

[8], but requires knowledge of the kurtosis of the stress response. Formulations for the 
kurtosis of the stress response in terms of the kurtosis of the load and the FRF were given for 
both the steady and the burst non-Gaussian load types. These were validated using time 
domain simulations of the stress response. Predictions of fatigue life using this new approach 
will be considered in future work. 

Finally, the excitations considered in this work were such that maximum stress range was 
below the yield stress of 400 MPa. However, higher loading levels could result in material 
and/or large-deflection nonlinearities. An assessment of the applicability of the spectral 
fatigue analysis, or development of some approximate form thereof, would need to be 
undertaken under such conditions. However, such an endeavour is beyond the scope of the 
present effort. 

ACKNOWLEDGMENTS 

The authors would like to acknowledge Paul Hollingsworth from Jaguar Land Rover, UK for 
the FE model. 

REFERENCES 

[1] Rizzi, S.A., Przekop, A., and Turner, T.L., "On the response of a nonlinear structure to 
high kurtosis non-Gaussian random loadings," EURODYN 2011, 8th International 
Conference on Structural Dynamics, Paper 41, Leuven, Belgium, July 4-6, 2011. 

[2] Newland, D.E., An introduction to random vibrations and spectral analysis (2nd 
edition). New York, Longman Inc, 1984. 

[3] Smallwood, D.O., "Vibration with non-Gaussian noise," Journal of the Institute of 
Environmental Sciences and Technology, Vol. 52, No. 3, pp. 13-30, 2009. 

[4] Sweitzer, K.A., "Random vibration response statistics for fatigue analysis of nonlinear 
structures," Ph.D. Dissertation, Institute of Sound and Vibration Research, University 
of Southampton, UK, 2006. 

[5] Papoulis, A., Probability, random variables, and stochastic processes (3rd edition), 
McGraw-Hill, Inc., 1991. 

[6] S.R., W., "Nonlinear vibration models for extremes and fatigue," ASCE Journal of 
Engineering Mechanics, Vol. 114, No. 10, pp. 1772-1790, 1988. 

[7] Ochi, M.K. and Ahn, K., "Probability distribution applicable to non-Gaussian random 
processes," Probabilistic Engineering Mechanics , Vol. 9, pp. 255-264, 1994. 

[8] Benasciutti, D. and Tovo, R., "Cycle distribution and fatigue damage assessment in 
broad-band non-Gaussian random processes," Probabilistic Engineering Mechanics, 
Vol. 20, pp. 115-127,2005. 

[9] Steinwolf, A., "Two methods for random shaker testing with low kurtosis," Sound & 
Vibration, pp. 18-22, October, 2008. 

[10] "ABAQUS online documentation, ABAQUS analysis user's manual." Providence, RI: 
Dassault Systemes Simulia Corp., 2012. 

[11] "MIL-STD-810F: Test method standard for environmental engineering considerations 
and laboratory tests," United States Department of Defense, 2000. 


[12] Matsuishi, M. and Endo, T., "Fatigue of metals subjected to varying stress," Japan 
Society of Mechanical Engineers, Fukuoka, Japan, March, 1968. 

[13] Halfpenny, A., "Rainflow cycle counting and fatigue analysis from PSD," 
Proceedings of the ASTELAB Conference, Paris, France, September 25 - 27, 2007. 

[14] "nCode DesignFife theory guide," HBM-nCode Products, Rotherham, UK 2009. 

[15] Falanne, C., "Mechanical vibration & shock (volume 4)," Hermes Penton Ftd, Fondon 

2002 . 

[16] Rice, S.O., "Mathematical analysis of random noise," Bell System Technical Journal, 
Vol. 23, pp. 282-332, 1944. 

[17] Downing, S.D. and Socie, D.F., "Simple rainflow counting algorithms," International 
Journal of Fatigue, Vol. 4, No. 1, pp. 31-40, 1982. 

[18] Palmgren, A., "The service life of ball bearings," NASA Technical Translation TT F- 
13460 (Translation of "Die lebensdauer von kugellagern," Zeitschrift des Vereines 
Deutscher Ingenierure, Vol. 68, No. 14, pp. 339-341, 1924), 1971. 

[19] Miner, M.A., "Cumulative damage in fatigue," Trans. ASME, Journal of Applied 
Mechanics, Vol. 67, pp. A159-A164, 1945. 

[20] Facoume, J.-F., Amblard, P.-O., and Comon, P., Statistiques d’ordre superieur pour le 
traitement du signal, Masson, 1997. 

[21] Papoulis, A., "Narrow-band systems and Gaussianity," IEEE Transactions on 
Information Theory, Vol. 18, No. 1, pp. 20-27, 1972. 

[22] Bendat, J.S. and Piersol, A.G., Random data: Analysis and measurement procedures, 
Wiley-Interscience, 1971. 

[23] Gardner, W.A., Introduction to random processes with applications to signals and 
systems (2nd edition), McGraw-Hill, 1989. 

[24] Roychowdhury, J., Fong, D., and Feldman, P., "Cyclostationary noise analysis of large 
RF Circuits with multi-tone excitation," IEEE Journal of Solid-State Circuits, Vol. 33, 
pp. 324-336, 1998. 


